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Abstract —  This  paper  describes  the  estimation  of  propagation 
loss  and  it  statistical  properties  based  on  observations  of  radar  sea 
clutter  data.  This  problem  is  solved  by  first  finding  an  ensemble 
(about  105  models)  of  relevant  refractivity  model  parameters 
and  then  using  all  these  models  to  map  into  the  propagation 
loss  domain.  In  this  mapping  each  refractivity  model  is  weighted 
according  to  its  data  likelihood  function. 

I.  INTRODUCTION 

Historically,  the  vertical  profile  of  atmospheric  refractivity 
of  a  “standard  atmosphere”  has  been  used  as  a  basis  for 
radar  performance  predictions.  Since  the  mid-1980’s  it  has 
become  common  practice  to  improve  this  profile  with  surface 
observations  of  air  temperature,  sea  temperature,  wind  speed 
and  relative  humidity,  in  conjunction  with  an  atmospheric 
surface-layer  model  to  estimate  the  refractivity  profile  within 
the  marine  atmospheric  surface  layer  [1].  The  refractivity 
profile  above  the  surface  layer  is  not  characterized  by  the 
surface  layer  models  and  is  usually  characterized  by  weather- 
balloons  or  rocket-sondes  [2].  But  whereas  the  surface  layer 
measurements  can  be  made  continuously,  weather-balloons  or 
rocket  soundings  are  expensive  in  both  use  of  expendables  and 
manpower.  An  alternative  means  for  estimating  refractivity  for 
the  purpose  of  radar  performance  assessment  is  to  use  radar 
clutter  (e.g.,  Fig.  1)  as  described  by  the  authors  in  [3],  [4]. 

A  general  weakness  of  all  methods  employed  in  estimating 
radar  performance  has  been  the  lack  of  a  means  for  quantifying 
the  impact  of  uncertainty  in  the  estimates  of  refractivity  that 
the  performance  estimates  are  based  on.  This  might  be  a 
particular  concern  when  using  radar  clutter  for  the  charac¬ 
terization.  Consider  Fig.  2:  A  pattern  of  radar  clutter  over  tens 
of  kilometers  on  a  given  azimuth  is  observed.  Unknown  to 
us,  however,  are  how  random  processes  in  individual  range 
bins,  and  colored  processes  in  range  (i.e.,  horizontal  gradients 
in  wind  speed,  hence  the  sea  clutter  radar  cross  section)  have 
affected  the  pattern.  In  the  context  of  our  forward  modeling, 
though,  this  makes  the  observation  somewhat  uncertain.  Many 
different  realizations  of  refractivity — possibly  from  disjoint 
regions  in  the  environmental  parameter  space — map  into  the 
region  of  uncertainty  surrounding  our  observation. 

The  principle  of  the  inversion  is  indicated  in  Fig.  2.  From  the 
observed  radar  clutter  data  or  the  data-space  d,  a  prediction  of 
propagation  loss,  the  information  usage  domain  u,  is  searched. 
The  vector  d  represents  the  obseved  radar  data  at  N  range 


Reflectivity  image:  April  02, 1998  Map  #  040298-12  18:00:00.3 
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Fig.  1.  Clutter  map  from  SPANDAR  corresponding  to  Wallops  Run  12. 


bins  and  the  vector  u  represents  a  vector  of  propagation  loss 
at  certain  height  and  range  bins.  This  involves  a  number  of 
steps  as  outlined  below: 

(1)  Determine  a  model  for  the  refractive  environment  (Ap¬ 
pendix  of  [3])  and  select  an  appropriate  propagation  model 
(e.g.,  TPEM  [5]).  These  two  form  the  mapping  d(m)  from  the 
parameter  space  m  to  data  space  d.  The  vector  m  represents 
the  unknown  environmental  parameters. 

(2)  Determine  the  mapping  u(m)  from  the  parameter  space 
to  information  usage  space.  Except  for  a  change  in  geometry 
(receiver  height),  this  is  similar  to  d(m). 

(3)  Find  acceptable  models  m  from  the  data.  As  indicated  in 
Fig.  2,  a  region  around  the  data  can  map  into  several  acceptable 
solutions  in  the  model  domain. 

(5)  Map  the  acceptable  models  into  the  data  usage  domain. 

(6)  The  posteriori  results,  both  in  parameter  m  and  data 
utility  u  domains  should  be  interpreted  probabilistically. 

As  indicated  in  Fig.  2,  this  mapping  is  non-unique.  There  are 
many  solutions  that  give  about  the  same  goodness-of-fit.  The 
maximum  likelihood  (ML)  is  the  best  fit  of  all.  In  particular 
we  can  associate  cost  (or  risk)  from  errors  in  the  solution — a 
reasonable  assumption  in  this  problem — and  the  uncertainty  of 
the  solution  is  non-uniform  in  the  data  usage  space.  Knowing 
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parameters  m  and  the  propagation  loss  usage  domain  u(m). 

the  posteriori  distribution  or  part  thereof  is  preferable  to  having 
a  single  point  estimate  such  as  the  ML  solution. 

A.  Replica  generation 

The  replica  clutter  return  d(m)  is  modeled  as  follows: 

1.  An  environmental  model  maps  the  environmental  parame¬ 
ters  m  into  modified  refractivity  as  a  function  of  range  r  and 
height  z. 

2.  A  combined  radar  and  propagation  model  maps  modified 
refractivity  into  radar  observations  d0(m). 

For  solving  of  the  inverse  problem,  a  large  number  (about 
90.000)  realizations  of  the  replica  have  to  be  generated. 

II.  Inverse  problem  framework 

In  the  Bayesian  paradigm,  the  solution  to  determining  pa¬ 
rameters  of  interest  m  given  an  observation  d  is  characterized 
by  the  posteriori  probability  pp(m  |  d).  A  clear  discussion  on 
on  inverse  theory  from  a  probabilistic  point  of  view  is  given  by 
[6] .  Additional  details  of  Monte  Carlo  sampling  of  posteriori 
probabilities  is  given  by  [7],  [8],  [9],  [10],  [11]. 

Following  Mosegaard  and  Tarantola  [10],  represent  the 
inverse  problem  as  conjuncture  of  information.  The  posteriori 
distribution  pp(m,d)  is  given  by 

P p(m,  d)  =  pe(m,  d)p0(m,  d)  ,  (1) 

where  p0(m,d)  is  the  prior  probability  and  pe(m,  d)  is  the 
probability  obtained  from  the  experiment.  We  note  that  this  is 
quite  similar  to  the  Bayesian  approach. 

The  joint  prior  theoretical  probability  is  constructed  as 
follows:  First  the  marginal  prior  probability  p(m)  is  used 
to  generate  samples  of  m  and  then  via  the  forward  model 
m)  the  joint  prior  probability  distribution  p(m,  d(m)).  It 
should  be  noted  that  this  is  not  just  a  simple  mapping  as  d(m) 
contains  additional  noise  and  shifting  of  the  clutter  in  range 
has  been  added  to  the  clutter  return  via  Markov  chains.  Note, 
that  the  construction  does  not  involve  any  observed  data  and 
thus  can  be  precomputed. 


The  experimental  probability  is  based  on  the  likelihood 
formulation  in  Sect.  III.  The  objective  function  0(d,  d(m)) 
is  a  squared  error  (Sect.  III). 

From  the  joint  posteriori  distribution  pp(m,d),  we  obtain 
the  posteriori  distribution  pp(m)  of  the  environmental  param¬ 
eters 

Pp(m)  =  [  P P(m,d)  dd  .  (2) 

Jv 

This  distribution  contains  all  relevant  information  and  from 
this  distribution  all  relevant  features  of  the  environment  can 
be  found  as  the  most  common  value  of  m  in  pp(m)  is  the 
maximum  a  posteriori  (MAP)  estimator. 

We  are  not  interested  in  the  environment  itself  but  rather 
better  estimates  of  the  parameters  in  the  information  usage 
domain  u(m).  Based  on  the  posteriori  distribution  pp(m),  the 
distribution  of  pp(u(m))  is  obtained  and  from  this  distribution 
all  relevant  statistics  of  the  usage  space  can  be  obtained. 

A.  Probability  of  utility  information  u 

Both  the  experimental  data  d  and  the  utility  data  u  are 
related  to  m  via  forward  models  F{ m)  and  G(m)  cor¬ 
responding  to  d(m)  and  u(m)  .  Thus  formally  we  have 
u  =  G(F~1( d)).  However,  this  direct  mapping  is  ill  posed 
and  is  instead  interpreted  based  on  conjuncture  of  information 
whereby  we  can  also  include  prior  information.  Similar  to 
Eq.  1  the  posteriori  probability  is  obtained 

Pp(u,m,d)  =pe(u,m,d)p0(u,m,d)  .  (3) 

Often  the  prior  probabilities  are  independent  and  the  above 
equation  can  be  simplified 

Pp(u,m,d)  =pe(u,m,d)pa(u)pa(m)pa(d)  .  (4) 

We  obtain  the  marginal  posteriori  distribution  pp(u) 

p  (u)  =  /  /  p  (u,m,d)  dmdd  .  (5) 

Jv  Jm 

Often  the  posteriori  distribution  of  pp(m)  is  first  obtained 
and  from  this  the  probability  distribution  of  the  propagation 
loss  is  obtained 

P P[u(r,z)}  =  j  u(r,z,  m)pp(m)  dm  (6) 

From  this  distribution  all  relevant  statistical  features  of  the 
propagation  loss  can  be  obtained,  as  for  example  the  mean 
and  variation  of  the  propagation  loss 

E  [u  ( r ,  z)\  =  J  u(r,  z)pp[u(r,  z)\  du  (7) 
Var  [u  ( r ,  z)]  =  J  (u  ( r ,  z)  -  E  [u(r,  z )])2  p p[u(r,  z)\  du  (8) 
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III.  Likelihood  and  Objective  function 

The  likelihood  and  the  objective  function  are  derived  as¬ 
suming  a  simple  linear  model 

d  =  Sd(m)  +  e  .  (9) 


Model  parameter 

Lower  bound 

Upper  bound 

Base  height,  z-^  (m) 

0 

400 

M-deficit,  M d  (M-units) 

0 

250 

Thickness,  ^thick  (m) 

0 

150 

Slope  (M-units/m) 

-1 

0.13 

Evaporation  duct  height  (m) 

0 

50 

the  scalar  5  represents  unknown  radar  strength  and  calibration 
constants,  d  and  d(m)  represents  measuered  and  modelled 
data  organized  into  a  vector  of  N  data  points,  corresponding  to 
each  range  bin.  m)  repensent  the  M  environmental  parameter 
to  be  determined.  Assuming  the  errors  e  to  be  Gaussian 
distributed  with  mean  Sd(m)  and  standard  deviation  C.  The 
errors  represent  all  features  that  are  not  modeled  in  the  data 
as  noise,  theoretical  errors,  and  modeling  errors.  This  is  a 
reasonable  model  for  high  signal-to-error  ratio,  but  for  low 
signal-to-error  ratio  this  is  not  a  good  assumption  [8].  The 
likelihood  function  is 

L(m)  =  (7r|C|)-iV  exp  (— [d  —  5d(m)]’i'C_1[d  —  Sd(m)])  . 

(10) 

where  N  is  the  number  of  data  points.  As  usual,  we  assume 
C  =  v\  The  source  5  can  be  estimated  in  closed  form  by 
requiring  Slog L/dS  =  0,  whereby 


5ml  — 


dM(m) 

||d(m||2 


function  is  then 

L( m)  =  {/kv)~n  exp  (— 0(m)/z/) 


TABLE  I 

Inversion  model  with  parameter  search  bounds.  These 

PARAMETER  BOUNDS  WERE  USED  BOTH  AT  THE  INITIAL  RANGE  (r  —  0) 
AND  AT  THE  FINAL  RANGE  ( r  —  lOOKM).  FOR  ALL  RANGES  EACH 
PARAMETER  WAS  LINEARLY  INTERPOLATED. 


p100 


where 


<f)(  m)  =  d2  — 


/dM(m)\2 

VlldMlJ  ’ 


(13) 


is  the  objective  function.  The  likelihood  estimate  of  the  noise 
z/ML  can  be  estimated  in  closed  form  by  solving  d  log  L/dv  = 


0, 


vML  =  4>(m)/N  . 


(14) 


Reinserting  this  into  the  likelihood  function  gives 


L(m) 


N 

7re0(m) 


N 


(15) 


The  ML  solution  mML  is  obtained  by  maximizing  the 
objective  function  over  all  m.  Finally,  an  overall  estimate  for 
the  error  power  v  is  obtained  from  (14)  at  the  environmental 
ML  solution:  z/ML(mML)  and  can  be  re-inserted  into  the 
likelihood  function.  For  simplicity,  we  consider  the  error  as 
known  and  only  keep  the  free  argument  m  of  the  objective 
function  0.  This  approach  leads  to  [8] 


£(m) 


N 


7T0(mML) 


N 


exp 


0(mML) 


(16) 


However,  posteriori  densities  also  can  be  based  on  (15). 
The  expressions  (15)  and  (16)  differ  in  their  error  estimates:  In 
(15),  the  error  power  is  estimated  for  each  geoacoustic  model 
vector  m.  In  (16)  the  global  ML  estimate  of  the  error  power 


Fig.  3.  The  refractivity  profile  (a)  used  to  generate  the  data,  and  (b)  the 
propagation  loss  data. 


(which  was  found  from  optimization)  is  posteriori  estimated 
and  then  is  (as  a  prior)  applied  to  all  model  vectors  m  [8]. 

The  above  derivation  assumes  that  the  error  in  each  sample 
is  uncorrelated  with  the  next  sample.  In  practice  these  are 
strongly  correlated  and  therefore  the  number  of  samples  N 
in  the  above  equations  must  be  replaced  with  the  effective 
number  of  samples,  A7eff . 

IV.  Example 

The  data  are  generated  based  on  the  helicopter  measured 
range-dependent  refractivity  profiles  (Run  7)  for  the  Wallops 
98  experiment  [3].  These  profiles  were  interpolated  using  the 
LARRI  program  [12].  A  range  interval  from  10-100  km  is 
used.  A  simple  trilinear  model  is  used  for  the  refractivity 
profile  as  outlined  in  the  Appendix  of  [3].  We  then  search  for 
refractivity  parameters  at  0  and  100  km  range  with  unknown 
parameters  as  given  in  Table  I.  To  obtain  refractivity  profiles  at 
other  ranges  the  parameters  are  interpolated  linearly.  The  first 
3  parameters  were  given  a  uniform  distribution  but  the  slope 
was  given  a  non-uniform  distribution  as  indicated  in  Fig  4 
(left  bottom).  This  is  because  a  negative  slope  is  only  likely 
for  shallow  ducts. 

Using  this  environmental  model,  90,000  replicas  are  pre¬ 
computed.  To  compute  the  prior  distribution  of  the  envi- 
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Fig.  4.  Prior  (left)  and  posteriori  (right)  distributions  of  the  parameter 
estimates. 


ronmental  parameters  all  samples  are  weighted  evenly  (left 
column  Fig.  4).  The  smoothness  of  these  indicate  that  the  prior 
is  well  sampled.  To  compute  the  posteriori  distribution  (right 
column  Fig.  4)  Eqs.  (3)  and  (16)  are  used.  We  note  that  the 
overall  behavior  seems  reasonable.  A  more  robust  convergence 
criterium  [9]  will  indicate  that  the  distributions  are  not  yet 
completely  converged. 

The  important  issue  is  to  be  able  to  estimate  statistics  of 
posteriori  propagation  loss.  Figs.  5-9  illustrate  this  estima¬ 
tion.  We  compute  the  priori  propagation  loss  based  on  an 
even  weighting  of  the  propagation  loss  from  each  generated 
refractivity  model.  The  posteriori  probability  distribution  of 
the  propagation  loss  is  based  on  weighting  the  propagation  loss 
from  each  refractivity  model  with  the  posteriori  probability, 
Eq.  (6). 

The  average  prior  propagation  loss  and  the  posteriori 
propagation  loss  are  shown  in  Fig.  5b  and  c.  It  is  seen  that  the 
posteriori  propagation  loss  identifies  a  ducting  environment  as 
observed  in  the  data,  Fig.  5a,  but  the  prior  does  not. 

The  probability  distribution  of  the  field  then  is  computed. 
Both  prior  field  (Fig.  6)  and  the  posterior  field  (Fig.  7)  are 
calculated  at  10  and  100  m  height  as  a  function  of  range. 
Figure  8  shows  the  probability  distribution  at  50  km  range 
and  10  m  height.  This  corresponds  to  a  simple  cut  though  the 
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Fig.  5.  The  propagation-loss  field  based  on  (a)  the  true  environment 
(from  Fig.  3a),  based  on  (b)  prior  information,  and  (c)  posteriori  information 
(bottom). 
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Fig.  6.  Prior  probability  distribution  for  propagation  loss  versus  range  at  10 
m  (top)  and  100  m  (bottom)  height. 


contour  plots  in  Figs.  6  and  7. 

Based  on  the  distribution  of  propagation  loss,  the  mean  and 
standard  deviation  are  computed  as  a  function  of  range  (10- 
100  km)  in  the  duct  (12-m  height)  or  above  the  duct  (102-m 
height),  see  Fig.  9.  Clearly,  the  posteriori  distribution  is  much 
more  compact  for  the  lower  height  than  for  the  higher  height 
—  a  behavior  that  is  expected  given  that  the  input  data  is 
effectively  a  (somewhat)  contaminated  direct  measure  of  the 
propagation  loss  at  the  surface. 

V.  SUMMARY 

An  algorithm  for  estimating  statistical  properties  of  prop¬ 
agation  loss  based  on  refractivity  from  radar  clutter  has 
been  described  using  a  likelihood  formulation.  The  likelihood 
function  is  developed  assuming  the  error  in  the  data  to  be 
Gaussian. 
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Range  (km) 

Fig.  7.  Posteriori  probability  distribution  for  propagation  loss  versus  range 
at  10  m  (top)  and  100  m  (bottom)  height. 


Range:  54.8045km 


0.25 


Propagation  loss  (dB) 

Fig.  8.  Posteriori  (solid)  and  prior  (dashed)  probability  at  50  km  range  and 
10  m  height. 
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